{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "4b61fc87",
   "metadata": {},
   "source": [
    "# 哈密顿量的构造\n",
    "*Copyright (c) 2021 Institute for Quantum Computing, Baidu Inc. All Rights Reserved.*"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "52532fe2",
   "metadata": {},
   "source": [
    "## 概览\n",
    "\n",
    "在本教程中，我们将展示如何使用 Paddle Quantum 的 `qchem` 模块通过化学分子式来构造适用于量子计算机的哈密顿量。我们将一步一步地学习如何从一个分子结构构造二次量子化的哈密顿量，以及如何将它转换成一组泡利矩阵。  \n",
    "\n",
    "哈密顿量是一个与物理系统总能量有关的物理量。一般来说，它可以表示为  \n",
    "\n",
    "$$\n",
    "\\hat{H}=\\hat{T}+\\hat{V},\\tag{1}\n",
    "$$\n",
    "\n",
    "其中 $\\hat{T}$ 代表动能，$\\hat{V}$ 代表势能。哈密顿量在变分量子算法中有着重要应用，如比[变分量子本征求解器](./VQE_CN.ipynb) 以及 [利用 Product Formula 模拟时间演化](./HamiltonianSimulation_CN.ipynb)。\n",
    "\n",
    "当我们用量子力学来解决化学问题时，需要写出一个哈密顿量来描述这个问题涉及的化学系统。我们可以根据这个哈密顿量计算该系统的基态和激发态，并利用这些信息进一步探索量子系统的所有物理性质。涉及电子结构问题的哈密顿量可以写成如下形式：\n",
    "\n",
    "$$\n",
    "\\hat{H}=\\sum_{i=1}^N\\left(-\\frac{1}{2}\\nabla_{x_i}^2\\right)+\\sum_{i=1}^N\\sum_{j< i}\\frac{1}{|x_i-x_j|}-\\sum_{i=1}^N\\sum_{I=1}^M\\frac{Z_I}{|x_i-R_I|},\\tag{2}\n",
    "$$\n",
    "\n",
    "该公式使用了[原子单位](https://en.wikipedia.org/wiki/Hartree_atomic_units)。该公式包含了 $N$ 个电子和 $M$ 个原子核，其中 $x_i$ 是第 $i$ 个电子的位置，$R_I$ 是第 $I$ 个原子核的位置。\n",
    "\n",
    "本教程分为四部分，首先我们先来讨论如何使用 `qchem` 模块构造一个分子。之后，我们将简要描述如何通过调用外部量子化学来计算\n",
    "[Hartree Fock](https://en.wikipedia.org/wiki/Hartree%E2%80%93Fock_method)\n",
    "单粒子轨道。接下来，我们展示如何得到二次量子化的哈密顿量。最后，我们将描述如何将费米子的哈密顿量 (Fermionic Hamiltonian) 转化为适用于量子计算机的泡利字符串 (Pauli strings)。  "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3160994c",
   "metadata": {},
   "source": [
    "## 水分子的结构\n",
    "\n",
    "在该部分，我们将展示如何从水的分子式和原子坐标构造出水分子。\n",
    "\n",
    "![h2o.png](figures/buildingmolecule-fig-h2o.png)\n",
    "\n",
    "在量桨中，我们将分子写成一个列表的形式，第一个元素是原子符号，第二个元素是它的位置 (Cartesian coordinate)。因此，用于表述分子的列表是由原子列表组成的。\n",
    "\n",
    "**注意：关于环境设置，请参考 [README_CN.md](https://github.com/PaddlePaddle/Quantum/blob/master/README_CN.md).**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "e6787ece",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Eliminate noisy python warnings\n",
    "import warnings\n",
    "\n",
    "warnings.filterwarnings(\"ignore\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "94b2032c",
   "metadata": {},
   "outputs": [],
   "source": [
    "# in Angstrom\n",
    "h2o_structure_direct = [[\"H\", [-0.02111417,0.8350417,1.47688078]],  # H 代表着水分子中的氢元素\n",
    "                        [\"O\", [0.0, 0.0, 0.0]],                     # O 代表着水分子中的氧元素\n",
    "                        [\"H\", [-0.00201087,0.45191737,-0.27300254]]]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "d354c5dd",
   "metadata": {},
   "outputs": [],
   "source": [
    "from paddle_quantum.qchem import geometry\n",
    "\n",
    "h2o_structure_xyz = geometry(file=\"h2o.xyz\")\n",
    "assert h2o_structure_xyz == h2o_structure_direct"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e996795a",
   "metadata": {},
   "source": [
    "## Hartree Fock 轨道的计算\n",
    "\n",
    "Hartree Fock 方法使用 [Slater 行列式](https://en.wikipedia.org/wiki/Slater_determinant)来表示 $N$ 电子波函数。它可以为我们提供一套单粒子轨道，这些轨道通常被作为其他量子化学方法的输入。 \n",
    "\n",
    "量桨使用 psi4 [1] 作为它的量子化学引擎。我们可以使用 `qchem` 模块提供的 `get_molecular_data` 函数来执行相关计算，得到分子的所需信息。 `get_molecular_data` 函数以分子结构、分子总电荷和自旋多重性为主要输入，并返回一个 OpenFermion [2] `MolecularData` 对象。  \n",
    "\n",
    "我们继续以水分子为例。为了运行 Hartree Fock 计算，我们需要将 `method` 设置为 *scf* (Self - Consistent Field)。 我们还可以通过在 `basis` 参数中指定 [basis set](https://en.wikipedia.org/wiki/Basis_set_(chemistry)) 的类型来提高 Hartree Fock 计算的精度。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "29945676",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Hartree-Fock energy for H2-O1_sto-3g_singlet (10 electrons) is -73.96770387867429.\n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.qchem import get_molecular_data\n",
    "\n",
    "h2o_moledata = get_molecular_data(\n",
    "    h2o_structure_direct,\n",
    "    charge=0,                # 水分子是中性的，不带电\n",
    "    multiplicity=1,          # 水分子只有一个未配对电子\n",
    "    basis=\"sto-3g\",\n",
    "    method=\"scf\",\n",
    "    if_save=True,            # 是否将 MolecularData 中的信息存储成 hdf5 文件\n",
    "    if_print=True,           # 是否需要打印出水分子的基态能量\n",
    "    name=\"\",                 # 指定 hdf5 文件的名字\n",
    "    file_path=\".\"            # 指定 hdf5 文件的路径          \n",
    ")\n",
    "\n",
    "from openfermion.chem import MolecularData\n",
    "\n",
    "assert isinstance(h2o_moledata, MolecularData)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "973af7d1",
   "metadata": {},
   "source": [
    "## 哈密顿量的二次量子化形式\n",
    "\n",
    "当我们研究电子系统时，通常会将哈密顿量写成[二次量子化](https://en.wikipedia.org/wiki/Second_quantization)的形式\n",
    "\n",
    "$$\n",
    "\\hat{H}=\\sum_{p,q}h_{pq}\\hat{c}^{\\dagger}_p\\hat{c}_q+\\frac{1}{2}\\sum_{p,q,r,s}v_{pqrs}\\hat{c}^{\\dagger}_p\\hat{c}^{\\dagger}_q\\hat{c}_r\\hat{c}_s,\\tag{3}$$\n",
    "\n",
    "其中 $p$, $q$, $r$ 和 $s$ 是 Hartree Fock 轨道，$\\hat{c}^{\\dagger}_p$ 和 $\\hat{c}_q$ 分别是生成算子 (creation operation) 和湮灭算子 (annihilation operation)。系数 $h_{pq}$ 和 $v_{pqrs}$ 分别是单体积分 (one-body integrations) 和双体积分 (two-body integrations)，该信息可以用如下的方式从 `MolecularData` 提取出来。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "a9a906ad",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-3.2911e+01  5.5623e-01  2.8755e-01  9.7627e-16 -7.4568e-02 -9.4552e-02  2.8670e-01]\n",
      " [ 5.5623e-01 -8.0729e+00 -4.0904e-02  4.2578e-16  1.7890e-01  3.5048e-01 -1.3460e+00]\n",
      " [ 2.8755e-01 -4.0904e-02 -7.3355e+00  1.1465e-16  4.1911e-01  5.2109e-01  7.0928e-01]\n",
      " [ 9.7627e-16  4.2578e-16  1.1465e-16 -7.5108e+00  4.1730e-15  8.3317e-15 -8.4993e-16]\n",
      " [-7.4568e-02  1.7890e-01  4.1911e-01  4.1730e-15 -5.7849e+00  2.0887e+00  1.2427e-01]\n",
      " [-9.4552e-02  3.5048e-01  5.2109e-01  8.3317e-15  2.0887e+00 -5.0803e+00  1.3967e-02]\n",
      " [ 2.8670e-01 -1.3460e+00  7.0928e-01 -8.4993e-16  1.2427e-01  1.3967e-02 -5.0218e+00]]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np \n",
    "np.set_printoptions(precision=4, linewidth=150)\n",
    "\n",
    "hpq, vpqrs = h2o_moledata.get_integrals()\n",
    "assert np.shape(hpq)==(7, 7)             # When use sto3g basis, the total number of molecular orbitals used in water calculation is 7\n",
    "assert np.shape(vpqrs)==(7, 7, 7, 7)\n",
    "\n",
    "print(hpq)\n",
    "# print(vpqrs)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "de4f1bf2",
   "metadata": {},
   "source": [
    "大多数情况下，我们并不需要把这些积分信息手动的提取出来，*qchem* 模块已经将该步骤融入函数 `fermionic_hamiltonian` 中，我们可以直接进行下一步的运算，获得哈密顿量。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "069cfcd5",
   "metadata": {},
   "outputs": [],
   "source": [
    "from paddle_quantum.qchem import fermionic_hamiltonian\n",
    "\n",
    "H_of_water = fermionic_hamiltonian(\n",
    "    h2o_moledata,\n",
    "    multiplicity=1,\n",
    "    active_electrons=4,\n",
    "    active_orbitals=4\n",
    ")\n",
    "\n",
    "from openfermion.ops import FermionOperator\n",
    "\n",
    "assert isinstance(H_of_water, FermionOperator)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "47f24b4b",
   "metadata": {},
   "source": [
    "通过指定 `active_electrons` 和 `active_orbitals`，我们可以限制哈密顿的自由度，从而可以减少哈密顿量的的项数，加速运算。我们还可以使用 *qchem* 中 `active_space` 函数来查看哪些轨道是**核心 (core)** 轨道，哪些轨道是**活跃 (active)** 轨道。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "c136344d",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "List of core orbitals: [0, 1, 2]\n",
      "List of active orbitals: [3, 4, 5, 6]\n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.qchem import active_space\n",
    "\n",
    "core_orbits_list, act_orbits_list = active_space(\n",
    "    10,                        # number of electrons in water molecule\n",
    "    7,                         # number of molecular orbitals in water molecule\n",
    "    active_electrons=4,\n",
    "    active_orbitals=4\n",
    ")\n",
    "\n",
    "print(\"List of core orbitals: {:}\".format(core_orbits_list))\n",
    "print(\"List of active orbitals: {:}\".format(act_orbits_list))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "05b6aefc",
   "metadata": {},
   "source": [
    "## 从费米子哈密顿量到自旋哈密顿量\n",
    "\n",
    "在量子计算中，我们只能使用由泡利矩阵构成的算符(operator)\n",
    "\n",
    "$$\n",
    "\\boldsymbol{\\sigma}_x=\\begin{pmatrix}\n",
    "0 & 1\\\\\n",
    "1 & 0\n",
    "\\end{pmatrix},\\quad \\boldsymbol{\\sigma}_y=\\begin{pmatrix}\n",
    "0 & -i\\\\\n",
    "i & 0\n",
    "\\end{pmatrix},\\quad \\boldsymbol{\\sigma}_z=\\begin{pmatrix}\n",
    "1 & 0\\\\\n",
    "0 & -1\n",
    "\\end{pmatrix}.\\tag{4}\n",
    "$$\n",
    "\n",
    "因此，我们需要将现有的哈密顿量（二次量子化形式）变换成量子比特算符 (qubit operator)，这种变换被称为 [Jordan-Wigner 变换](https://en.wikipedia.org/wiki/Jordan%E2%80%93Wigner_transformation)。\n",
    "\n",
    "> 此外，Bravyi-Kitaev 变换也可以达到相同的效果，只需要将参数 mapping_method 改成 'bravyi_kitaev' 即可。\n",
    "\n",
    "在量桨中，我们提供 `spin_hamiltonian` 函数来实现从费米子哈密顿量到自旋哈密顿量的变换，代码如下："
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "7a831c80",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "There are 193 terms in H2O Hamiltonian in total.\n",
      "The first 10 terms are \n",
      " -72.10615980544183 I\n",
      "-0.007310917992546845 X0, X1, Y2, Y3\n",
      "0.005246087073083395 X0, X1, Y2, Z3, Z4, Y5\n",
      "0.0016283548447088131 X0, X1, Y2, Z3, Z4, Z5, Z6, Y7\n",
      "0.005246087073083395 X0, X1, X3, X4\n",
      "0.0016283548447088131 X0, X1, X3, Z4, Z5, X6\n",
      "-0.005994544380559041 X0, X1, Y4, Y5\n",
      "0.0013876441781026563 X0, X1, Y4, Z5, Z6, Y7\n",
      "0.001387644178102656 X0, X1, X5, X6\n",
      "-0.009538223793221256 X0, X1, Y6, Y7\n"
     ]
    }
   ],
   "source": [
    "from paddle_quantum.qchem import spin_hamiltonian\n",
    "\n",
    "pauli_H_of_water_ = spin_hamiltonian(\n",
    "    h2o_moledata,\n",
    "    multiplicity=1,\n",
    "    active_electrons=4,\n",
    "    active_orbitals=4,\n",
    "    mapping_method='jordan_wigner'\n",
    ")\n",
    "\n",
    "print('There are', pauli_H_of_water_.n_terms, 'terms in H2O Hamiltonian in total.')\n",
    "print('The first 10 terms are \\n', pauli_H_of_water_[:10])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f5af604d",
   "metadata": {},
   "source": [
    "现在，我们知道了如何从分子结构构造出相应的哈密顿量，接下来不妨去看一看如何使用[变分量子本征求解器](./VQE_CN.ipynb) (VQE)去计算氢分子的基态能量。"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d7417802",
   "metadata": {},
   "source": [
    "---\n",
    "## 参考文献\n",
    "\n",
    "[1] [Psi4: an open-source ab initio electronic structure program](https://wires.onlinelibrary.wiley.com/doi/abs/10.1002/wcms.93)\n",
    "\n",
    "[2] [OpenFermion: the electronic structure package for quantum computers\n",
    "](https://iopscience.iop.org/article/10.1088/2058-9565/ab8ebc)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
